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D. I. Pontin*, A. K. Al-Hachami* & K. Galsgaard* 



Abstract 



o 

(N 

oo 
6 

o3 



> 
00 

\o 

in 
o 



X 



Context: Three-dimensional (3D) null points are present in abundance in the solar corona, and the same is likely to 
be true in other astrophysical environments. Recent results from solar observations and from simulations suggest 
that reconnection at such 3D nulls may play an important role in the coronal dynamics. 

Aims: The properties of the torsional spine and torsional fan modes of magnetic reconnection at 3D nulls are 
investigated. New analytical models are developed, which for the first time include a current layer that is spatially 
localised around the null, extending along either the spine or the fan of the null. The principal aim is to investigate 
the effect of varying the degree of asymmetry of the null point magnetic field on the resulting reconnection process 
- where previous studies always considered a non-generic radially symmetric null. 

Methods: Analytical solutions are derived for the steady kinematic equations, and are compared with the results of 
numerical simulations in which the full set of resistive MHD equations is solved. 

Results: The geometry of the current layers within which torsional spine and torsional fan reconnection occur 
is strongly dependent on the symmetry of the magnetic field. Torsional spine reconnection occurs in a narrow 
tube around the spine, with elliptical cross-section when the fan eigenvalues are different. The eccentricity of the 
ellipse increases as the degree of asymmetry increases, with the short axis of the ellipse being along the strong field 
direction. The spatiotemporal peak current, and the peak reconnection rate attained, are found not to depend 
strongly on the degree of asymmetry. For torsional fan reconnection, the reconnection occurs in a planar disk in 
the fan surface, which is again elliptical when the symmetry of the magnetic field is broken. The short axis of the 
ellipse is along the weak field direction, with the current being peaked in these weak field regions. The peak current 
and peak reconnection rate in this case are clearly dependent on the asymmetry, with the peak current increasing 
but the reconnection rate decreasing as the degree of asymmetry is increased. 



1 Introduction 

There has been significant progress in recent years towards understanding the properties of magnetic reconnection in 
three dimensions (3D). In particular, it is now appreciated that the intense current layers necessary for reconnection 
in astrophysical plasmas may form at a number of different characteristic structures of the magnetic field. One such 
structure is a 3D magnetic null point. Such null points have been demonstrated to be present in abundance in the solar 



corona (e.g. Longcope & Parnell 2009), and the same is highly likely to be true in other astrophysical environments 
such as the coronae of other stars and of accretion disks. Recent results from observations and simulations suggest 



that reconnection at these nulls may play an important role in the coronal dynamics (e.g. Luoni et al. 2007 Lynch 



et al. 2008 Pariat et al.l 2009 Masson et al. 2009). Note also that the importance of reconnection at 3D nulls is not 



restircted to astrophysical plasmas, but plays a role both in the Earth's magnetosphere (e.g. Xiao et al. . 2006 ) and 
some laboratory plasmas (e.g. Gray et al. 2010). 



Recent studies have revealed a number of characteristic modes of reconnection that may occur at such 3D nulls. 
These have recently been categorised by Priest & Pontin ( 2009[ ) into 'torsional spine reconnection', 'torsional fan 
reconnection' and 'spine-fan reconnection'. Torsional spine and torsional fan reconnection occur when an equilibrium 
magnetic null point field is disturbed by a rotational perturbation (the rotation being in a plane perpendicular to the 
spine). Spine- fan reconnection occurs when a shear perturbation is applied that disturbs the locations of the spine 
and fan - this leads to a localised current layer forming around the null and flux transport across the fan separatrix 
surface. 
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Close to the null the magnetic field may be written 



B = VB x 



(1) 



and the eigenvalues and eigenvectors of the matrix VB determine the locations of the spine and fan of the null (see e.g. 
Parnell et al. 1996). Previous studies of the generation of current layers at 3D nulls due to rotational motions have 



considered only the case where the background equilibrium null point has a rotational symmetry, i.e. in which the two 
eigenvalues associated with the fan are equal. During the resulting evolution in which torsional spine and torsional 
fan reconnection takes place, the respective current layers, plasma flows, etc., have exhibited azimuthal symmetry 
due to the azimuthally symmetric background magnetic fields and perturbations. We have previously investigated the 
effect of varying the symmetry of the background null point field on the spine-fan reconnection mode, and shown that 
while the topological properties of the reconnection are unchanged, the dimensions and intensity of the current layer, 



and the reconnection rate, are strongly dependent on the degree of asymmetry (Al-Hachami & Pontin 2010), as well 
as the relative angle between the shear driving and the strong/weak field directions (Galsgaard & Pontin 20111. In 
this paper we generalise existing models for torsional spine and torsional fan null point reconnection as follows. We 
begin in each case by introducing a new kinematic analytical solution for the corresponding reconnection mode in 
which a localised current layer is present at the null. We then go on to consider the effect of varying the symmetry 
of the background magnetic field, by varying the ratio of the fan eigenvalues of the null. Lastly we perform numerical 
simulations of the full system of MHD equations in which null point configurations with varying degrees of symmetry 
are subjected to rotational disturbances, to complement the analytical models. 

The paper is organised as follows. In Sect. [2] we review previous modelling of torsional spine and torsional fan 
reconnection. In Sect. [3] we introduce our new analytical and numerical models for torsional spine reconnection, and 
discuss their implications. We go on in Sect. [4] to do the same for torsional fan reconnection. Finally in Sect. [5] we 
present our conclusions. 



2 Existing models for torsional spine and fan reconnection 



Both torsional spine and torsional fan reconnection involve the formation of a current layer in which the current vector 
is directed parallel to the spine line at the null. It was first shown by Pontin et al. (2004) that the corresponding 



reconnection takes the form of a rotational slippage of magnetic flux threading the non-ideal region. (This is in contrast 
to the case where the current vector is parallel to the fan, in which case magnetic flux is reconnected across the spine 



and fan (Pontin et al. 2005).) The magnetic flux undergoes this rotational slippage in response to rotational flows in 



the ideal region in which the rate of flux transport in the azimuthal direction is different for field lines entering the 
non- ideal region than it is for field lines exiting the non- ideal region. The original model of Pontin et al. (2004) is 
based on the magnetic field 

B = B [r,jr/2,-2z] 

in cylindrical polar coordinates. This results in a spatially uniform current parallel to the spine (z-axis). In order 
to obtain a 3D-localised non-ideal region (relevant in astrophysical plasmas which are approximately ideal almost 
everywhere), it was therefore necessary to impose a localised resistivity profile. 

In the following sections we introduce two new analytical kinematic solutions. In these solutions we have been 
able, for the first time, to include a current density which is fully spatially localised with its peak intensity focussed at 
the spine or fan. This current localisation allows a spatially uniform resistivity to be used, adding a degree of physical 
plausibility to the models, since in practice in an astrophysical plasma a localised non-ideal region is associated with a 
localised current layer. In each of our solutions, the structure chosen for the magnetic field and the resulting current 
layer is based on the results of resistive MHD simulations in which the dynamic formation of these current layers has 
been observed. It is worth noting that additional kinematic solutions with a spatially varying current density have 
been presented by Wyper & Jain (2010[). In their solutions the current is spatially localised, but in contrast to our 



new solutions is focussed away from the null point. 

Resistive MHD simulations have demonstrated that the form of the current layer is different depending on whether 
the rotational perturbation primarily disturbs the fan field lines or field lines around the spine. The perturbation 
behaves essentially as an Alfven wave, travelling along the magnetic field lines, which due to the hyperbolic structure 
of the field around the null leads to the perturbation accumulating either in the vicinity of the spine or the fan. When 
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the fan field lines are subjected to a rotational disturbance, torsional spine reconnection takes place in a tubular current 



structure that forms at the spine (Rickard & Titov 1996 Pontin & Galsgaard 2007]). Within this tube, the magnetic 
field lines spiral around the spine line. The radius of the tube decreases, and the current intensifies, until the twisting 
of the field being driven by the perturbation is balanced by rotational slippage facilitated by resistive diffusion. The 
reconnection rate, defined in 3D as the maximal value of 

* = J E-B/|B|dZ (2) 

along any field line threading the non-ideal region, then quantifies this rotational slippage. 

When field lines in the vicinity of the spine line are disturbed, a current layer forms on the fan surface leading to 



torsional fan reconnection (Rickard & Titov 1996 Galsgaard et al. 20031. Again field lines spiral within the current 



layer, whose magnitude intensifies as the twisting of the field is concentrated in an increasingly narrow sheet over the 
fan surface. Once the sheet becomes sufficiently thin resistive diffusion dissipates the twist leading once again to a 
rotational slippage of field lines. 



3 Torsional spine reconnection 



3.1 Torsional spine reconnection: kinematic model 

In this section we describe a kinematic model for torsional spine reconnection within a current tube that is localised to 
the spine of the null point. The form of the magnetic field and resulting current structure is chosen to match behaviour 



observed in the numerical simulations described by Pontin & Galsgaard ( 2007 ) . The method of solution is the same as 



that described by Hornig & Priest (2003) and Pontin et al. (2004). Specifically, we solve the steady-state, kinematic, 



resistive MHD equations 



E + v x B = T]3 
V x E = 



VxB = fi J 
VB = 0. 



(3) 



We specify the form of the magnetic field B and solve for the electric field (E) and plasma flow perpendicular to B 
(v_l) via 

(E — 77J) x B 



77J • B ds, 



E 



V_L 



B 2 



(4) 



(the component of v parallel to B being arbitrary). Note that s is a parameter along magnetic field lines, i.e. the 
integration of $ is performed along these field lines. The resistivity rj = t/q is taken to be uniform. 



We consider first the rotationally symmetric case. We set 

Bo 



B 



L 



(Bp + Bj) 



where 



Bj = 



2// 



Bp = [r,0,-2z] 

[0,0,0] 



r < a 
& \z\ < b 

otherwise 



(5) 

(6) 
(7) 



in cylindrical polar coordinates (r, 9, z) where j, a, b are positive real numbers and n, m, fi, v are positive integers. Bp 
defines the potential 'background' null point component, and Bj defines the non-potential component of B associated 
with the tubular current structure, which extends to radius r = a and to z = ±6. In accordance with previous 
simulation results we assume an extended tube of current aligned to the spine so that 6>a. For m = n = fj, = v = l, 
B and J are continuous and differentiable. However, in order that all physical quantities in the final solution are 
continuous and differentiable we find that it is necessary to choose m — 3,/i = n = 2 and v = 1. The current density 
is represented by the plot in Fig. [l] 
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Figure 1: Left Magnetic field lines for the torsional spine model defined by Eqs. ( 5|8|9 ) for a = 1, b = 5, j = 
3,p = 1.5, q — 1. The shaded surface shows a current isosurface. Right, current vectors in the y = plane for 
a=l,b = 2,q = l. 



(|5]|7|) can be found in a straightforward way 
is independent of 9 we require only 



Parametric equations for the magnetic field lines associated with Eqs. 
by solving 9X/<9s = B(X(s)). Since the integrand in the first equation in (I4J 
r = r{) exp(_Bos/Lo) & n d z = zoexjp(—2Bos/Lo) (it is straightforward to also obtain an expression for 9(tq, 9q, Zq, s)). 
Solving Eqs. Q with s = at z = ±zo we find that the qualitative properties of the electric field and flow field are 
similar to those discovered by Pontin et al. (2004). In particular, rotational plasma flows are still present around the 
spine axis of the null. Using the freedom of arbitrary flow parallel to B in the model we can for illustrative purposes 
choose to add a component vn such that v z = 0. We then see a purely azimuthal flow. Due to the return current 
region close to r — a (see Fig. [T]) the sense of rotation of the plasma changes at an intermediate radius, as shown in 
the top-right frame of Fig. [2j 

We now investigate how the properties of the solution vary when the rotational symmetry of the above system is 
broken. When the rotational symmetry is lost it is no longer possible to find closed-form expressions for the field line 
mapping. We therefore numerically integrate B to find field lines and solve Eqs. ^ on a rectangular grid. We may 
break the symmetry either in the potential component Bp defining the magnetic null or in the component Bj defining 
the current tube. Our new potential component of the magnetic field is given by 



B f 



2p 
l X 'p+l 



y, 



-2z 



(8) 



in Cartesian coordinates where p > is a parameter. As p varies the magnetic field along the spine direction is fixed 
while the ratio between the fan eigenvalues (associated with the eigenvectors along the x and y directions) varies (see 
Fig. [lj. We choose to break the symmetry in Bj by converting to Cartesian coordinates and setting 



B 



/ — 



i(i-(f) 6 ) 4 (i-(f) 4 ) a [-w.«.o] 

[0,0,0] 



R < a & \z\ < b 
otherwise 



(9) 



where R 2 — x 2 +qy 2 (note that this reduces to expression Q when q = 1). This has the effect of distorting the current 
into a cylinder with elliptical cross-section, with major and minor axes along the x- and y-axes, extending to x = ±a, 
y = ±a/y/q (see the images on the left in Fig. [2j). Pre-empting the results of the following section, we present here 
results for p — q, such that as p increases the current tube narrows along the direction associated with the large fan 
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eigenvalue, i.e. the strong field direction in the fan. We set So = Lq = rjo = j = 1, a = 1, b = 4, and solve Eqs. Q on 
a rectangular grid with 81 gridpoints in each direction covering the volume —2 < x, y<2,0<z<4 with the solution 
being symmetric about z = 0. We restrict our attention to the range p > 1, which simply selects the y direction as 
the strong field direction in the fan. 

The results of the above analysis are presented in Figs. [2][4j As p is increased, the current tube shrinks in the 
y-direction, with the dominant current component J z intensifying in the part of the tube close to the y-axis (i.e. the 
direction of the short axis of the ellipse) - see Fig. [2] The stronger current in this region results in an enhanced plasma 
flow speed. The direction of the flow is also distorted from the circular pattern at p = 1, but continues to flow on 
closed elliptical paths around the spine (z-axis). As the fan plane is approached the radius of the elliptical shells of 
positive and negative azimuthal flow increase, owing to the hyperbolic nature of the field structure - see Fig. [3] 

In order to determine the reconnection rate we calculate VP as defined in Eq. ([2|. Due to the breaking of the 
symmetry it is no longer clear that the maximal value of should occur along the spine field line, as was found in 
previous studies (note that the current modulus has maximum value away from the spine for large p) . However, it turns 
out that indeed the maximum occurs along field lines asymptotically close to the spine for all p. Figure [3] displays 
the peak value of the current density (which we impose) and the reconnection rate as a function of the degree of 
asymmetry. It is clear that the peak current scales linearly with p(— q) and that correspondingly the reconnection rate 
scales linearly with p. Note however that all of the above solutions were obtained with a fixed value of the parameter j 
which also contributes to controlling the peak current density, and that the velocity and reconnection rate will increase 
proportional to this parameter. As we shall see below, however, we do not have a linear increase of the peak current 
with p in our simulations, and it may therefore be more realistic to take j = j(p). 



3.2 Torsional spine reconnection: simulations 



To complement the above kinematic model we now describe the results of numerical simulations of the full system 
of resistive MHD equations. The code has been extensively used for other similar simulations, the details of the 



scheme being described by Nordlund & Galsgaard (19971. The MHD equations are non-dimensionalised by setting 



the magnetic permeability /xq = 1 and the gas constant equal to the mean molecular weight. The result is that for a 
volume in which \p\ — |B| = 1, time units are such that an Alfven wave would travel one space unit in one unit of 
time. The equations are solved on a numerical grid of 256 3 gridpoints over [x, y, z] G [±0.75, ±0.75, ±3] with uniform 
i] = 10~ 4 . The grid is stretched to give higher resolution around the null, where Ax = Ay = 0.0034, Az = 0.014. All 



boundaries are closed and line-tied. We repeat the simulations described by Pontin & Galsgaard (2007), in which a 



localised rotational perturbation of the magnetic field is imposed on a background equilibrium null point. Specifically, 
we begin with a potential magnetic field given by Eq. (|8|), plasma velocity zero everywhere, and initialise the plasma 
density and thermal energy to be spatially uniform with values p = 1 and e = 0.025, respectively. As a result there 
exist regions of both high and low plasma-/3 in the domain. (3 approaches infinity as one approaches the null, and 
takes a value of 1 on an ellipsoidal surface which cuts the coordinate axes at x = ±0.09(p +1), y = ±0.09(p + l)/p, 
z = ±0.09. 

In addition to the above, we impose at t = a magnetic field perturbation composed of a ring of magnetic flux 
centred on the null point and lying in the fan plane: 



B e = b exp 



(r - r ) 



C 2 



(10) 



with b = 0.05, r = 0.18, a = 0.08, ( = 0.06 (see Figs. 1 and 2 of |Pontin & Galsgaard| ( |2007[ )). (We have also performed 
simulations where this initial perturbation is elliptical rather than circular, but found no change to the qualitative 
results - thus here we confine our discussion to the circular perturbation.) 

For t > the perturbation splits, with wavefronts travelling both toward and away from the null. We focus 
on the behaviour of the ingoing pulse, which for p — 1 gradually stretches out to form a cylindrical tube of intense 
current around the spine. When p ^ 1, the azimuthal symmetry of the perturbation wavefront is broken as soon as 
the evolution begins. The Alfven speed in the radial direction now depends on x and y, and so the wavefront travels 
toward the spine faster along the y-axis where |B| is stronger. Propagation in the z-direction is essentially unaffected, 
and the current distribution associated with the perturbation forms into a cylinder with elliptical cross-section, whose 
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Figure 2: Current density (left) and plasma flow vectors (right) for the torsional spine kinematic solution. Plotted in 
the plane z = 2 for Bo = Lq = j = a = 1,6 = 4 and p = q = 1 (above), p = q = 2 (middle) p = q = 5 (below). The 
grey curve (red online) in the right-hand plots marks the boundary of the non-ideal region. 
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Figure 3: Intensity map of v x in the x = plane for the kinematic torsional spine solution with p = q = 1. 




Figure 4: Dependence on the anisotropy parameter q = p of the maximum values of |J| and the reconnection rate 
^max for the kinematic torsional spine model. 



Table 1: Data on the simulations of torsional spine reconnection. 



p 


I J Imarc 




V]/ a 

1 max 




Lx b 


T b 


L z b 


1 


0.31 


0.24 


4.2 x 10" 


5 


0.03 


0.03 


3.9 


2 


0.29 


0.27 


5.7 x 10" 


5 


0.09 


0.023 


3.9 


3 


0.31 


0.31 


6.0 x 10" 


5 


0.13 


0.019 


3.3 


5 


0.34 


0.34 


6.1 x 10" 


5 


0.19 


0.018 


2.5 


10 


0.36 


0.36 


5.8 x 10" 


5 


0.26 


0.017 


1.9 



peak integrated parallel electric field attained. 



b current layer dimensions measured at the time when J z > reaches its temporal maximum. 
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t=0 



t=0.9 



t=1.8 



t=3.0 



Figure 5: Time sequence showing contours of J z in the z = (left) and x = (right) planes, for the torsional spine 
simulation with p = 2. 

length and eccentricity both increase as the pulse steepens towards the spine and the peak current correspondingly 
intensifies (see Figs. [5j [6]). 

As shown in Fig. [7j there is a strong (at least 100%) increase in the peak current density as the disturbance reaches 
the spine. It is also clear from the plot that this occurs earlier for larger p due to the increased speed of propagation 
along the y-axis. In general, the peak current is higher for larger p. However, this is complicated by a competing effect 
- namely that there are two distinct maxima in |J| max during the simulations, one corresponding to the localisation 
of the leading edge of the pulse and the other to the trailing edge of the pulse. As can be readily seen in Fig. [5j they 
correspond to opposite signs of J z . Thus, while the maximum positive value of J z strictly increases as p increases, 
it is found that for p — 1,2 the trailing edge of the pulse (corresponding to J z < 0) dominates. The dimensions of 
the current layer are shown in the final three columns of Table [l] They are measured at a time corresponding to the 
localisation of the leading edge of the pulse, i.e. the time when J z > reaches its temporal maximum (which occurs 
at a time marked by the peak on each curve in Fig. [7] between t = 1 and t = 2, except for the run with p = 1 in which 
J z reaches its maximum at t = 2.1). The dimensions of the current tube in the rry-plane demonstrate its elliptical 
nature (centred on the spine) with eccentricity increasing with p as shown in Fig. [6j The length in the z-direction 
(L z ) decreases with increasing p. This is because the z-component of Bp (and thus the propagation speed in z) is 
independent of p, and because the pulse localises in the xy-plane at an earlier time for larger p as explained above. It 
is also worth noting from Fig. [6] that for p ^ 1 the maximum current is attained not exactly on the spine, but in two 
locations displaced symmetrically from the spine along the £-axis. 

We turn now to consider the reconnection rate, calculated as described in Eq. The maximal value of ^ is 
found over all field lines that thread the current layer, these being field lines that pass close to the null and its spine 
and fan. For p = 1, due to symmetry, the maximum can be found on any field line which runs asymptotically close to 
the null. However, for p > 1 the current is strongest in the weak field regions around the x-axis, and so field lines that 
pass through these regions attain the highest values of The evolution of the reconnection rate is plotted for the 
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-0.2 -0.1 0.0 0.1 0.2 0.3 



Figure 6: Above: Current isosurface at 50% of maximum value, for the torsional spine simulations with p = 1 (left), 
p = 2 (middle) and p = 5 (right). Taken in each case at the time when the positive value of J z reaches a maximum. 
Below: contours of J z in the z = (fan) plane over [x, y] £ [±0.2, ±0.2] for the same runs at the same times. 



0.40 




Figure 7: Evolution of the peak value of |J| for the different torsional spine simulations: p = 1 (solid line), p = 2 
(dotted), p = 3 (dashed), p = 5 (triple-dot-dashed) and p = 10 (dot-dashed). 
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Figure 8: Evolution of the reconnection rate for the different torsional spine simulations: p — 1 (solid line), p = 2 
(dotted), p = 3 (dashed), p = 5 (triple-dot-dashed) and p — 10 (dot-dashed). 

different simulations in Fig. § The maximum value attained does not depend strongly on p (see the fourth column 
of Table [T]), except that it is significantly lower for p = 1. However, this is likely down to the fact that for p = 1 the 
disturbance extends all the way to the z-boundaries, and so we 'miss' some of the length of the current layer. Note 
that the reconnection rate calculated here is a net effect of integrating through regions of both positive and negative 
E • B - we have different senses of reconnection (rotational slippage) occurring on the leading and trailing edges of the 
pulse, and here we measure the net effect. 



3.3 Torsional spine reconnection: discussion 

The numerical simulations discussed above demonstrate the localisation of a rotational perturbation towards the spine 
of a non-symmetric linear 3D null. The resulting current intensification is associated with torsional spine magnetic 
reconnection. The current tube that forms around the spine has a structure that is closely matched by the kinematic 
steady-state model described in Sect. |3.1| In particular, the current is dominated by the component parallel to the 
spine (J z ), and is localised within a tube of elliptical cross-section. The short axis of the ellipse is aligned with the 
strong field direction in the fan plane, along which the current is most intense. The eccentricity of the ellipse increases 
as the magnetic field asymmetry increases. This dynamic evolution should be compared with the stationary analytical 
solution at the stage where the current layer has become fully localised on the spine, as in the analytical solution 
- i.e. when the maximum (or minimum) of J z is achieved. If it were feasible to drive the perturbation from the 
boundaries in a steady fashion, we would expect this spine-localised current tube to persist in a quasi-steady state for 
a period controlled by the period of steady driving, as in Sect. [4] below - see Fig. 13 (driving a rotational motion in a 



smooth way around the fan boundaries in our Cartesian domain is not practical). 

The plasma flow also has a similar qualitative structure in the kinematic model and simulations. This structure is 
that of a rotation along elliptical paths around the spine, with these elliptical paths closely following the current density 
contours. One difference between the model and simulations is that both signs of rotational flow are only seen during a 
certain period in the simulations. Specifically, as the twist associated with the perturbation propagates towards the null 
it drives flow predominantly in the positive rotational sense, and when the reconnection process is completed (t > 3.5) 
the field then 'untwists' leading to a large-scale flow in the opposite direction. It is only approximately between the 
times that J z reaches its maximum positive and negative values (approximately 1.5 < t < 3, see Fig. [7| that the 
Lorentz force accelerates the plasma in opposite rotational senses at different distances from the null. By contrast, 
in the steady model where momentum balance is neglected, rotational plasma flows of both senses are required to 
maintain the steady state. 

In the numerical simulations it is found that the peak reconnection rate is approximately independent of p. We 
can understand why this should be the case by noting that, while J z increases with increasing p, the length of the 
current layer along z decreases with increasing p. Thus for larger p we have a larger integrand in Eq. ([2]), but it is 
non-zero over a shorter distance. Note also that the reconnection rate is found to be given by the integral of Eh along 
field lines lying asymptotically close to the spine. Since J is dominated by J z and the current layer has minimal extent 
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in the xy-p\ane, then it is natural that the reconnection rate should not depend strongly on the B Ky field components 
away from the null. In order to match these simulation results for the p-dependence of the reconnection rate, the 
magnetic field Bj in our kinematic model could be normalised by a factor proportional the current modulus. 

It is worth considering here the interpretation of the reconnection rate as measured above. Within the framework of 
general magnetic reconnection (Schindler et al. 19881, the flux reconnected is measured by calculating the differential 



transport of flux on opposite sides of the non-ideal region, for a flux surface passing through the peak of the pseudo- 
potential ^ (see also Hesse et al. 2005). This theory is based on the assumption that no nulls are present within 
the non-ideal region, while here of course we consider a non-ideal region containing a null point. In the rotationally 
symmetric case the maximal value of |^| is found on any field line passing asymptotically close to the null. This 
measures the differential rate of flux transport through any flux surface bounded by the spine and fan in the half-space 
(above or below the fan), by symmetry (Pontin et al. 2004 Pontin 2011). That is, it measures the net magnetic 



flux reconnected at the given time. However, in the present simulations for p ^ 1 the value of \I/ no longer exhibits 
azimuthal symmetry. Therefore the maximal value of |^| obtained above measures the maximum differential rate of 
flux transport for any flux surface bounded by the spine and fan, though one should bear in mind that this now varies 
around the spine. Of course in a steady-state situation, this differential rate of flux transport must be the same for all 
such flux surfaces (by continuity), and indeed it is found to be in the steady-state kinematic solution. 

The behaviour of the disturbance - spreading out along the spine in all of the simulations - suggests that the 
dominant wave mode is an Alfven wave. However, it is highly likely that other wave modes are present. Owing to 
the differing waves speeds approaching the null for different directions, one may speculate that if the disturbance 



were sufficiently strong as to be considered non-linear, then some of the effects discussed by McLaughlin et al. (2009) 
would be present. Those authors considered the propagation of a non- linear fast-mode wave towards a 2D X-point, 
and observed the formation of cusp-shaped structures as the wavefront collapsed, and the formation of both fast and 
slow mode shock waves. However, since our simulation is three-dimensional, it is not feasible at present to use the 
numerical resolution required to properly resolve such features. Note further that the perturbation in our simulations 
is initialised in the vicinity of the /3 — 1 surface. It is known that wave-mode conversion may take place in the vicinity 
of this surface. Thus the propagation of the disturbance towards the null, as well as the final form of the current layer, 
may be different if we were to repeat the series of simulations with different values of the plasma-/?. 



4 Torsional fan reconnection 



4.1 Torsional fan reconnection: kinematic model 



We now turn our attention to modelling of the torsional fan reconnection mode, which involves rotational slippage 
of field lines in a current layer localised around the fan surface. We proceed to solve Eqs. (13) in the same way as 
described in Sect. [3~Tj Again, we first analyse a model for the cylindrically symmetric case, in which for the first time 
a localised current layer is included in the fan plane. The structure of the magnetic field is chosen by comparing with 
the numerical simulations of Galsgaard et al. (2003); Pontin & Galsgaard ( 2007[ ). We again construct our magnetic 
field as the sum of a potential part (Bp) and non-potential part (Bj) as in Eq. ([5|, with Bp = [r, 0, —2z], and this 



time 



B, 



2 m 



2/i 



i 2)1 



r < a 
& \z\ < b 

otherwise, 



(11) 



[0,0,0] 

see Fig. [9] Note that Bg is now odd in z, and since we are modelling a current layer focused on the fan plane we 
assume that 6<Co. We again choose the integers m, n, [x, v in such a way that all physical quantities in our solution 
are continuous and differentiable, specifically m — 3, ll = n = 2, v = 6. The field line mapping for B is again given 
by r = ro exp(£?os/Lo) and z = zq exp(— 2Bqs/Lq), along with a lengthy expression for 6(ro, 9o, Zq, s) which is not 
required to obtain the solution. This time we solve Eqs. Q setting s — at r — tq — a. We choose to do this 
because it leads to non-zero flow for \z\ > which is consistent with the observed result from the simulations that 
this reconnection mode is set up by rotational driving flows in the regions around the spine footpoints. 

The new torsional fan solution is represented in Figs. |10| |11| The plasma flow has a rotational structure, as in 
the solution of Pontin et al. (2004). That is, when we subtract a component of v parallel to B such that v z = 0, 
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y 



Figure 9: Above: Magnetic field lines for the torsional fan model defined by Eqs. ( 5|8|12 1 for a = 5,b = 1 



50, p = 1.5, q — 1. The shaded surface shows a current isosurface. Below: current vectors in the y = plane for 
a = 3, b = 1, q = 1. 



then the remaining flow is non-zero only in the azimuthal direction (see the top-right frame of Fig. 10 1. Thus, field 
lines traced from comoving footpoints in the ideal region at z > z (or z < —z ) rotate around the spine at a fixed 
radius, while field lines traced from the ideal region at r > tq remain fixed (v = there), and we have a change of 
connectivity in the form of a counter-rotational slippage. Owing to the fact that J r changes sign for different levels of 
z, the rotational flow within the current layer has regions of both clockwise and anti-clockwise rotation, in a similar 
way to the torsional spine model discussed above. 

As before, we now investigate the dependence of the properties of our solution on the symmetry of the magnetic 
field. We proceed as in Sect. 3.1 to solve Eqs. Q using the semi-analytical method described there. 81 3 gridpoints 



are used over the domain —4 < x, y < 4, < z < 2, and we use parameters Bq = Lq = j '■ — 770 = 1, b = 1, a — 4. 

We take the potential component of our magnetic field (Bp) to be given by Eq. with the non-potential 
component taken to be 

^(i-(f) e ) 4 ( 1 -tt) 1 )' 2 [-»-' I -°] t1,f s % (12) 

[0, 0, 0] otherwise 



where R 2 = qx 2 + y 2 (which reduces to Eq. (Ill when q = 1). The current layer now has the shape of an elliptical 
disc, with major and minor axes along the x- and y-axes, extending to x = ±a/y/q, y = ±a. Pre-empting the results 
of the following section, we present here results for p = q, such that as p increases the current disc shrinks along the 
direction associated with the small fan eigenvalue, i.e. the weak field direction in the fan. As shown in Fig. [TUJ the 
current density is enhanced in the regions around the short axis of the ellipse. As in the torsional spine solution, if 
we set the parallel flow in such a way as to eliminate v z , then the plasma flow in the xy-plane follows closed elliptical 
paths, being strongest in magnitude where the current is enhanced. 

Examining the dependence of the peak current and reconnection rate on the degree of asymmetry, we find that 
both increase with increasing q, as shown in Fig. |11| 

4.2 Torsional fan reconnection: simulations 



We now perform numerical simulations similar to those described in Sect. 3.2 designed to investigate the effect of 
the background field asymmetry on the torsional fan reconnection mode (see Table [5| . The setup of the simulations 
is the same as before, except that this time we take 256 3 gridpoints over [x, y, z] € [±2.5, ±2.5, ±0.5]. This time we 
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Figure 10: Current density (left) and plasma flow vectors (right) in the kinematic torsional fan model. Plotted for 
Bq = L = j ' = b = 1, a = 4 and p = q = 1 (above), p — q = 2 (middle) p = q = 5 (below). The grey curve (red 
online) in the right-hand plots marks the boundary of the non-ideal region. The current density is plotted at z = 0, 
while the plasma flow is plotted in the plane z = 0.5. 
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Figure 11: Dependence on the anisotropy parameter q — p of the maximum values of |J| and the reconncction rate 
^max f° r the kinematic torsional fan model. 



Tab le 2: Data on the simulations of torsional fan reconncction. 



P \Jxy\max ^ 'max 



1 


0.94 


12.7 x 10 


-5 


2 


1.20 


7.4 x 10" 


5 


3 


1.33 


5.7 x 10" 


5 


5 


1.49 


4.5 x 10- 


5 


10 


1.55 


3.8 x 10- 


5 



perturb the initial equilibrium (Eq. ([8 



by applying a rotational driving velocity on the boundaries around the spine 

The sense of rotation is opposite above and 



footpoints both above and below the null, as in Galsgaard et al. (20031 
below the fan plane. Specifically, we apply the following azimuthal velocity in the z = ±0.5 planes 



Vg = A 



1 r(l + tanh(l - CV)) t < 2r 



(13) 



where r — yjx 2 + y 2 , C = 10, r = 1.6 and A = =p0.1 at z = ±0.5. We perturb the system using a different method to 



Sect. 3.2 for the following reasons. First, the rotation of the fan plane (Sect. 3.2 1 is done with an internal perturbation 
within the domain since driving from the boundary in this way in our Cartesian simulation domain would involve 
driving around the corners of the domain which is difficult to do in a smooth and stable way. Second, we choose to 
drive from the boundary in this section since applying an internal perturbation it is difficult to develop a significant 
current increase at the fan. This is because owing to the spatially uniform resistivity (which we cannot significantly 
reduce due to the limitations of numerical resolution) an internal perturbation (which must be well localised initially) 
suffers significant diffusion before reaching the fan, as discussed by|Pontin fc Galsgaard (2007). Note that this is not 



such an issue for the case where the fan field lines are rotated (Sect. 3.2) since in that case the incoming pulse spreads 
only in one direction (z) owing to the hyperbolic nature of the field while it contracts in two (x, y) - however, for the 
case of a spine rotation the in-coming pulse spreads in two directions {x,y) and contracts only in one (z). 

For the symmetric case p = 1, the evolution of the system following the initiation of the driving velocity is described 
in detail by Galsgaard et al. (2003). A disturbance that is dominated by a torsional Alfven wave is launched from 



the driving z boundaries towards the null. The wave front spreads along the x- and y-directions travelling along the 
field lines, with its velocity in the z-direction being independent of x and y. The wavefront steepens as it approaches 
the fan surface due to the hyperbolic structure of the field, eventually forming a planar current layer in the fan. This 



process is demonstrated in Fig. 12 (see also the left-hand image in Fig. 14), from which one can see that the current 



associated with the current layer is oriented parallel to the fan surface and flows radially inwards toward the null (note 
that in this x = plane, J x is zero by symmetry). It is worth noting that there also remains a large distributed current 
near the boundaries where the twisting was applied - observed as strong concentrations of J z - with the modulus of 
the current in this region being approximately equal to that of the current in the layer at the fan. We have repeated 
the simulation with different values of rj, and found that for lower r\ the peak current density at the fan is stronger 
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Figure 13: Evolution of the maximum of \J X y\ f° r the different torsional fan simulation runs: p — 1 (solid line), p = 2 
(dashed), p = 3 (dotted), p = 5 (dot-dashed) and p = 10 (triple-dot-dashed). For comparison, the time variation of 
the amplitude of the driving velocity is also represented, with the crosses. 




0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 

Ml 



Figure 14: Contours of |J| for the torsional fan simulations, plotted in the z = (fan) plane over [x, y] € [±2.5, ±2.5] 
for p = 1 (left), p — 2 (middle) and p = 5 (right). Taken in each case at the time when the value of \J X y\ reaches a 
maximum. 



relative to the concentration near the boundary. Hence, if we were able to run the simulation with a realistic value of 
r] for an astrophysical plasma the current layer at the fan would dominate. 

As a result of the dominance of the J z component near the driving boundaries (as discussed above), it is most 
clear to observe the formation and evolution of the fan current layer by plotting the evolution of the peak value of 



|Ja;j/|, as in Fig. 13 (the dominant current component in the current layer is parallel to the fan surface for all p). We 
can see that for p = 1 this quantity rises steadily as the pulse approaches the fan, after which there is a period of 
around t = 1.5 when the value remains steady, after which it decreases. Note that this steady period is consistent 
with the period during which the driving velocity remains steady - see the crosses in the figure. It is clear that the 
overall qualitative behaviour is similar between the simulations. However, the overall peak current attained is largest 
when the initial field is most asymmetric (p — 10). Note that for the simulations with the largest values of p the peak 
current does not remain steady for such a long period as for small p. This may be influenced by the fact that the 
disturbance interacts with the y-boundaries at later time for large p. 

We now examine the spatial structure of the current distribution for simulation runs with different p. When the 
boundary driving is initiated the azimuthal symmetry of the wavefront that propagates into the domain is broken. 
The wavefront remains approximately planar (z constant), spreading into an elliptical shape with long axis along the 
y-direction (along which the Alfven speed is greater) and short axis along x. However, as the wavefront gets closer to 
the null, it becomes steadily more inhomogeneous, with the current density focussing in the weak field regions near 
the x-axis. That is, although the current layer that forms eventually at the fan is more extended along y (as a simple 
consequence of mapping a circular driving region on the boundary along B), this current is most intense along the 
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Figure 15: Evolution of the reconnection rate for the different torsional fan simulation runs: p = 1 (solid line), p = 2 
(dashed), p = 3 (dotted), p = 5 (dot-dashed) and p = 10 (triple-dot-dashed). 



short axis of the ellipse, as shown in Fig. 14 As p is increased, the peak current increases as described above, and this 



increasing current maximum is localised in a gradually narrower 'channel' around the a;-axis. 

We turn now to consider the rate of reconnection in the different simulations. As discussed above, at the limited 
magnetic Reynolds number we are able to use the fan current layer does not dominate over the distributed current near 
the driving boundaries - though indications are that it would for more realistic astrophysical parameters. Therefore, 
in order to examine only the effect of reconnection in the thin current layer, we calculate the reconnection rate by 
integrating E\\ along the magnetic field line in the fan plane that passes through the peak of the current density - 
which in practice lies very close to the x-axis by symmetry. When the current layer has fully formed at the fan this 
genuinely measures the reconnection rate associated with the dynamically forming fan current layer. Clearly at earlier 
stages when the perturbation is yet to reach the fan there will still be a modest amount of reconnection, which is not 



measured by this approach. Therefore, when analysing the plots of reconnection rate versus time shown in Fig. 15 
one should bear in mind that the curves do not accurately portray the reconnection rate prior to t w 2.5. We see 
that the peak reconnection rate associated with the fan current layer occurs in each of the simulations at t 4, just 
before the peak current density starts its decline. A clear pattern emerges that the reconnection rate is greatest for 
the symmetric case with p — 1, and steadily decreases for simulations with larger p. 



4.3 Torsional fan reconnection: discussion 

The numerical simulations discussed above demonstrate the localisation of a rotational perturbation towards the fan 
of a non-symmetric linear 3D null, resulting in torsional fan magnetic reconnection. The planar current structure that 
forms around the fan qualitatively resembles that found in the kinematic steady-state model described in Sect. |3.1| 
In particular, the current is dominated by a component parallel to the fan where the current vector is directed (at 
its maximum intensity when p ^ 1) radially towards the null. Current isosurfaces have an elliptical shape, with the 
short axis of the ellipse aligned with the weak field direction in the fan plane, along which the current is most intense. 
The eccentricity of the ellipse increases as the magnetic field asymmetry increases. The plasma flow also has a similar 
qualitative structure in the kinematic model and simulations. This structure is that of a rotation along elliptical paths 
around the spine, with these elliptical paths closely following the current density contours. 

In the torsional fan reconnection simulations - in contrast to the torsional spine case - the peak reconnection 
rate depends strongly on the magnetic field asymmetry (p parameter). This is because in the torsional fan case, the 
contribution to the integrand in the reconnection rate definition ^ comes largely from field lines closely aligned to the 
xy-plane, which is the plane in which the magnetic field varies as p is varied. The interpretation of the reconnection 
rate when p ^ 1 should be the same as that described above in Sect. |3.3[ Perhaps counter-intuitively, the scaling of 
the reconnection rate with p is opposite to the dependence of the peak current on p (as p increases the peak current 



increases while the reconnection rate decreases - see Figs. 13 15 1. This is because, unlike in 2D, the reconnection rate 



is not a local quantity defined at a point, but rather is defined as an integral along a field line. So if the peak value 
of the integrand in Eq. ^ is larger, but the integrand is non-zero over a shorter section of the field line, this may 
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result in a lower net value for the integral, or vice-versa. Note that this behaviour was also observed in the study of 



spine-fan reconnection carried out by Al-Hachami & Pontin (20101. Note also that the effect was not captured in the 
kinematic solution - however, it is worthwhile noting that the scaling of the reconnection rate with p is much weaker 
than that of the peak current. Specifically, while the current increases by a factor 3 between p = 1 and p = 10 in 
the kinematic solution, the reconnection rate increases only by a factor w 1.2. 

In the simulations after the current layer that formed dynamically at the fan has started to dissipate, the distributed 
current near the driving boundaries dominates to a steadily greater extent, since the large-scale current concentrations 



dissipate only very slowly. This is consistent with the results of Pontin & Craig ( 2005 1 , who found that only when 



shear perturbations of the spine/fan are made does the lowest energy state of the system involve a current sheet at 
the null - when rotational perturbations are performed around the spine then the lowest energy state is achieved 
when the twist of the field is distributed along the field lines, in concentrations that extend outwards from the driving 
boundaries. 

Finally, note that here we have chosen to model an instance of torsional fan reconnection in which the vorticity of 
the driving flow around the spine has opposite sign on the opposite boundaries. As a consequence, a strong current 
- directed predominantly in the radial direction - is present at the fan surface. However, in our kinematic model we 
could equally have chosen Bg to be even in z to model the situation where the driving flows have the same sign of 



vorticity, which leads to cancellation at the fan plane of the currents generated (Galsgaard et al. 2003). It should also 
be noted that the return currents present close to z = ±b are not present in the simulations as we drive only in one 
sense, but could be induced by reversing the sign of the driving velocity at some intermediate time in the simulations. 



5 Conclusion 

Here we have presented analytical and numerical models for torsional spine and torsional fan magnetic reconnection at 
3D null points. The analytical models included for the first time fully localised current layers - focussed at the spine 
or fan - that determine the boundary of the non-ideal region, thus alleviating the requirement in previous models to 
have an artificially localised ('anomalous') resistivity. We also for the first time investigated the generic case where 
the null point is not radially symmetric, i.e. where the fan eigenvalues are not equal. 

3D null points have been demonstrated to be present in abundance in the solar corona, and the same is likely 
to be true in other astrophysical environments. Recent results from solar observations and from simulations suggest 
that reconnection at such 3D nulls may play an important role in the coronal dynamics. The fan separatrix surfaces 
of these 3D nulls divide the coronal magnetic field into distinct topological domains, i.e. distinct regions of magnetic 
connectivity between the photospheric flux concentrations. The torsional spine and torsional fan reconnection modes 
do not act to transfer magnetic flux between these topological domains, but rather permit a rotational slippage within 
the domains of the magnetic flux lying close to the nulls and therefore the domain boundaries. These null point 
reconnection modes are therefore unlikely to be involved in energetic events that involve a large-scale restructuring of 
the magnetic flux between topological domains as the coronal field seeks a lower energy state (e.g. during solar flares). 
Rather, they are a mechanism to dissipate energy and reduce stress associated with the dynamic perturbation of the 
coronal field by the turbulent boundary driving from the photosphere. 

We have shown that the geometry of the current layers within which torsional spine and torsional fan reconnection 
occur is strongly dependent on the symmetry of the magnetic field defining the null point. Torsional spine reconnection 
still occurs in a narrow tube around the spine, but with elliptical cross-section when the fan eigenvalues are different. 
The eccentricity of the ellipse increases as the degree of asymmetry increases, with the short axis of the ellipse being 
along the strong field direction. Furthermore, the current profile is not azimuthally symmetric around the spine, but is 
peaked in these strong field regions. The numerical simulations suggest that the spatiotemporal peak current, and the 
peak reconnection rate attained, do not depend strongly on the degree of asymmetry. For torsional fan reconnection, 
the reconnection occurs in a planar disk in the fan surface, which is again elliptical when the symmetry of the magnetic 
field is broken. The short axis of the ellipse is along the weak field direction, with the current being peaked in these 
weak field regions. The peak current and peak reconnection rate in this case are clearly dependent on the asymmetry, 
with the peak current increasing but the reconnection rate decreasing as the degree of asymmetry is increased. 

While we have relaxed the rotational symmetry of the magnetic field in these studies, the field structure of a 
linear null is still relatively simple. In the future it will be important to understand how these reconnection modes 
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are modified - and how they release the energy associated with imposed stresses - when the null point is embedded 
in a more realistic coronal geometry. One recent study suggests that other features present in the magnetic field may 
attract the current preferentially over the nulls and therefore may in some cases inhibit the formation of the torsional 
spine and torsional fan current layers (Santos et al. 2011 ). Other studies have focussed on explaining explosive events 



and thus have been motivated by identifying reconnection via flux transfer between topological domains, facilitated by 

In practice one would 



Masson et al. 


2009 


Pariat et al. 


2009 


Torok et al. 


2009) 



expect a combination of the pure reconnection modes to appear, and significant further study is required to determine 
the manifestations of, interactions between, and relative importance of the spine-fan mode and the two torsional 
modes discussed here. The importance of a number of other parameters such as the magnitude of the perturbation, 
the plasma-/3 and the resistivity rj is also yet to be explored. 
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